home *** CD-ROM | disk | FTP | other *** search
/ Aminet 33 / Aminet 33 - October 1999.iso / Aminet / misc / math / TCalcStats2c.lha / TCalcStats2c / AREXX / Deviates.rexx < prev    next >
Encoding:
OS/2 REXX Batch file  |  1999-07-03  |  3.9 KB  |  149 lines
  1. /* Calculate random Normal Deviates */
  2. options results
  3. if ~show('P','TCALC') then do
  4.     address command 'run turbocalc:turbocalc'
  5.     address command 'waitforport TCALC'
  6.     loadflag=1
  7. end
  8. address 'TCALC'
  9. 'DEFPUBSCREEN()'
  10. /* Add-in Rexx Libraries needed for some routines */
  11. signal on syntax
  12. if ~show('l','rexxmathlib.library') then
  13.    call addlib('rexxmathlib.library',0,-30)  /* add to library list */
  14. if ~show('l','rexxreqtools.library') then
  15.    call addlib('rexxreqtools.library',0,-30)
  16. if ~show('l','rexxsupport.library') then
  17.    call addlib('rexxsupport.library',0,-30)
  18. signal off syntax
  19.  
  20. /* Start Main Routine */
  21. if loadflag=1 then 'Load()'
  22. 'ActivateWindow()'
  23. NCols=rtgetlong("0","Enter number of columns to provide","Input Request") /*,,'rt_pubscrname="TCALC"')*/
  24. if rtresult=0|NCols=0 then do
  25.     'Message "Aborting!"'
  26.     'DEFPUBSCREEN("Workbench")'
  27.     exit
  28. end
  29. NRows=rtgetlong("0","Enter number of rows to provide","Input Request") /*,,'rt_pubscrname="TCALC"')*/
  30. if rtresult=0|NRows=0 then do
  31.     'Message "Aborting!"'
  32.     'DEFPUBSCREEN("Workbench")'
  33.     exit
  34. end
  35. /* Get cell reference for output range */
  36. out_cell=rtgetstring(,"Enter Cell Reference for Output","Input Request") /*,,'rt_pubscrname="TCALC"')*/
  37. if out_cell="" then do
  38.     'DEFPUBSCREEN("Workbench")'
  39.     exit
  40. end
  41. if length(out_cell)<2 | datatype(left(out_cell,1),'n') then do
  42.     'Message "Invalid cell reference"'
  43.     'DEFPUBSCREEN("Workbench")'
  44.     exit
  45. end
  46. /* Suppress Screen Redraw to Speed Things Up */
  47. 'Refresh 0'
  48.  
  49. /* Open a small output window on tcalc screen*/
  50. fo=0
  51. CR='0a'x
  52. DisplayMsg="Calculating...Please Wait."||CR||"User input is disabled during calculation."||CR
  53. if open(6Info, 'con:100/0/450/80/Progress/SCREEN TCALC', w) then do
  54.      call writeln(6Info, DisplayMsg)
  55.     fo=1
  56. end
  57. else do
  58.     'Message "TCALC Screen not available for Progress messages"'
  59. end
  60. CALL DELAY(150)
  61.  
  62. /* Create data */
  63. RN.=0
  64. Call RANDN(NCols,NRows)
  65. call writeln(6Info,"Writing output to window...")
  66.  
  67. /* Output */
  68. 'SelectCell' out_cell
  69. 'ColumnWidth 10'
  70. 'Put "Table of Random Normal Deviates with zero Mean and unit Variance"'
  71. 'CursorDown 1'
  72. 'GetCursorPos'
  73. top_cell=result
  74. Do x=1 to NCols
  75.     Do y=1 to NRows
  76.         'Put' RN.x.y
  77.         'CursorDown 1'
  78.     end
  79. 'SelectCell' top_cell
  80. 'Column' x
  81. end
  82. 'SelectCell' out_cell
  83. 'CursorDown' NRows+3
  84. 'Put "Calculated using the ratio of uniforms method"'
  85. 'CursorDown 1'
  86. 'Put "of A.J. Kinderman and J.F. Monahan"'
  87. 'CursorDown 1'
  88. 'Put "augmented with quadratic bounding curves."'
  89. 'CursorDown 1'
  90. 'Put "ADAPTED FROM WORK PUBLISHED IN TRANSACTIONS ON MATHEMATICAL SOFTWARE,"'
  91. 'CursorDown 1'
  92. 'Put "VOL. 18, NO. 4, DECEMBER, 1992, PP. 434-435."'
  93. 'Refresh 1'
  94. 'Refresh 2'
  95. DisplayMsg="Cleaning up ...."||CR||"Exiting"
  96. result=writeln(6Info, DisplayMsg)
  97. if result~=0 then do
  98.     /*Wait 3 seconds*/
  99.     CALL DELAY(150)
  100.     /* close window*/
  101.     result=close(6Info)
  102. end
  103. 'DEFPUBSCREEN("Workbench")'
  104. exit
  105.  
  106. /* Procedures */
  107.  
  108. RANDN: Procedure Expose RN.
  109.  
  110.     arg N,R
  111.     S=0.449871
  112.     T=-0.386595
  113.     A=0.19600
  114.     B=0.25472
  115.     R1=0.27597
  116.     R2=0.27846
  117.     Do i=1 to N
  118.       Do j=1 to R
  119.         Do forever
  120.             U = RANDU(time('s'))
  121.             V = RANDU(time('s'))
  122.             V = 1.7156 * (V - 0.5)
  123.         /*Evaluate the quadratic form*/
  124.             X  = U - S
  125.             Y  = ABS(V) - T
  126.             Q  = X**2 + Y*(A*Y - B*X)
  127.         /*Accept P if inside inner ellipse*/
  128.             IF (Q < R1) Then Leave
  129.         /*Reject P if outside outer ellipse*/
  130.             IF (Q > R2) Then Iterate
  131.         /*Reject P if outside acceptance region*/
  132.             IF (V**2 >  -4.0*LOG(U)*U**2) Then Iterate
  133.         End
  134.         RN.i.j=V/U
  135.       End
  136.     End
  137. Return
  138.  
  139. /*      ALGORITHM 712, COLLECTED ALGORITHMS FROM ACM.*/
  140. /*      THIS WORK PUBLISHED IN TRANSACTIONS ON MATHEMATICAL SOFTWARE,
  141.      VOL. 18, NO. 4, DECEMBER, 1992, PP. 434-435.
  142.   The function RANDN() returns a normally distributed pseudo-random
  143.   number with zero mean and unit variance.  Calls are made to a
  144.   function subprogram RANDU() which must return independent random
  145.   numbers uniform in the interval (0,1).
  146.  
  147.   The algorithm uses the ratio of uniforms method of A.J. Kinderman
  148.   and J.F. Monahan augmented with quadratic bounding curves.
  149. */